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We study the S=^ quantum spin ladder system with the four-spin exchange, using density 
matrix renormalization group method and an exact diagonalization method. Recently, the 
phase transition in this system and its universality class are studied. But there remain 
controversies whether the phase transition is second order type or the other type and the 
nature of order parameter. There are arguments that the massless phase appears. But this 
does not agree with our previous result. Analyzing DMRG data, we try a new approach in 
order to determine a phase which appears after the phase transition. We find that the edge 
state appears in the open boundary condition, investigating excitation energies of states with 
higher magnetizations. 

PACS numbers: 75.10.jm 

I. INTRODUCTION 

S=l/2 quantum spin two-leg ladder systems have been studied from both theoretical and 
experimental points of view [1]. These systems are related to high-T c superconductors and 
Haldane's conjecture [2]. For example, Dagotto et al. suggested that the superconductivity 
occurs in hole-doped ladder systems [3]. 

Usually, two spin exchanges have been discussed in quantum spin systems. By the way, 
spin exchange interactions originate from electron exchanges. For example, Heisenberg 
model is derived from the second order perturbation in the strong coupling limit of the 
Hubbard model. In contrast, higher order perturbations give many body spin exchange in- 
teractions, for example, (Sj ■ Sj) (S& ■ S;). These type interactions in two leg ladder systems 
are discussed by Nersesyan et al[4] . Investigating correlation functions, they referred to level 
k = 2 SU(2) Wess-Zumino-Witten (WZW) model or SU(2) c = 3/2 conformal field theory 
(CFT)). 

Recently four spin cyclic exchanges attract attention [5-14, 16, 18-20]. And from experi- 
ments, it has been suggested that a four spin cyclic exchange interaction plays an important 
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role in some systems, one dimensional two-leg ladder in La 6 Ca 8 Cu 2 404i [9, 10] and two 
dimensional square lattice in La 2 Cu0 4 [11]. The latter is known as a parent insulator of 
high-T c superconducting systems. The importance of a cyclic four spin exchange in such 
systems is pointed out by Honda et al [12]. In two-leg ladder systems (see Fig. 1), four 
spin exchanges play important roles. Sakai et al [13] and Nakasu et al [14] studied this 
system with magnetizations in a magnetic field, in respect of the magnetization plateau and 
plateau-gapless transition. This phase transition is the BKT type[15]. 

In our previous paper, we numerically studied s=l/2 spin ladder system under zero field, 
and we found a second order phase transition [16]. We studied the critical behavior of this 
system based on conformal field theory [17]. We calculated the central charge and one of 
the scaling dimension numerically, and we obtained ~ |. Therefore this phase 

transition is described as the c = § CFT or the k = 2 577(2) WZW model. This is a 
phase transition with Z 2 symmetry breaking. We discussed that the ordered state has the 
translational symmetry breaking and the rung-parity symmetry breaking, which supports 
the staggered dimer order (see Fig. 4). Recently this result is obtained by Miiller et al 
[18] in analytically. They used bosonization technique and its Majorana representation 
in the weak coupling limit. With the help of Pade- approximation in the strong coupling 
limit, they obtained the critical exponent for the energy gap, i] ~ 1 numerically. This 
corresponds with our result [16]. Lauchli et al [19] discussed a phase diagram for this model. 
They calculated staggered dimer structure factor with exact diagonalization and a local 
(S Q) j • S aj j+i) expectation value on one of the two leg in the dimer long range order phase 
with the density matrix renormalization group (DMRG) method. As a result, they conclude 
that there is a phase transition between a rung-singlet phase and a staggered dimer phase. It 
is clear and understandable to see the order parameter directly as Lauchli. But it is needed 
to inspect the accuracy of numerical calculation and finite size effects since the correlation 
length is large. In fact, as we can see later, we need a careful analysis on the DMRG data. 

However these results do not correspond to Honda et al [20]. They insisted that a massless 
phase appears in the strong four-spin coupling region. Honda et al investigated a spin 
gap and a spin-pair correlation function using DMRG method [21]. They found that at 
Jmng = Jieg — 1 an d J r ing — 0.3 a spin gap between the singlet (S^ ot = 0) and the triplet 
(^tot = 0' ^l) vanishes. 

Hikihara et al [22] discussed the self-duality and the chirality using Runge-Lenz vector. 
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And they referred to an incommensurate character. They derived the duality transformation 
which leaves the form of the Hamiltonian unchanged. This self-dual line is same to an exact 
line which was found by Kolezhuk et al [23]. With the matrix product method, one can 
exactly solve this model on this line [23]. The exact ground state is a product of rung- 
singlets on this line. 

Hikihara et al calculated spin gaps, a triplet excitation and a quintet excitation, using 
DMRG in open boundary condition. They found out that both gaps are quite small in J r i ng 
large region. They found that triplet gaps are exactly zero within numerical accuracy in 
some system sizes, but quintet gaps seem to be finite. They said that it is very difficult to 
estimate the critical point accurately due to the very slow vanishing of the spin gaps around 
the phase transition point. 

The extended massless phase is not consistent with some gapful phase. We will reanalyze 
the spin ladder system with four spin cyclic exchange. Miiller et al and Lauchli et al insisted 
that there is a phase transition from a gapped phase (rung-singlet phase) to another gapped 
phase (staggered dimer phase). Hikihara et al discussed that the quintet gap is quit small 
but it seems to be finite. Honda et al insist that there is a phase transition from a gapped 
phase (rung-singlet phase) to a gapless phase. But these results are not compatible. At least 
one interpretation is mistaken. 

We think that Honda et al found the edge state, since they used DMRG method with the 
open boundary condition. They only calculated the energy gap between the singlet state 
and the triplet state. So, we calculate energy gaps between the singlet and triplet, and 
between the singlet and the quintet. In this paper, we try to check the consistency between 
the second order phase transition and DMRG data using a new method. We try to calculate 
a correlation length quantitatively and to discuss the finite size effect. 

II. SPIN HAMILTONIAN AND FOUR SPIN CYCLIC EXCHANGE 

In this paper, we discuss the following Hamiltonian, 





(1) 



4 



where P, P^ 1 are four spin cyclic exchange operators (see Fig. 2), a = 1, 2 are indexes of spin 
chains (see Fig.l), J rung is an interchain coupling, Ji eg is an intrachain coupling, and J r i ng 
is a four spin cyclic exchange coupling. 

This four spin cyclic exchange operator can be expressed using spin operators as follows, 

Pi,i+i + Pi,i+i — ^ + Si^-S^j+i + S 2i j-S 2 ,i+i + S^j^^+i 

+ S2,j-Sl 5 i+l + Si,i-S 2 ,j + Si 5 j + i-S2,i+l 

+4{(Si ) j •S 2i j)(Si i j + i-S2.i+i) 

+ (Sl,i-Si,i+i) (S2,i-S 2i j + i) 

— (Si,i-S 2l i+i) (S 2 ,i-Si,i+i)}. (2) 

This interaction have been known from old times [24, 25]. 

Four spin coupling terms, such as (S a ,j • Spj) (S 7 ,jt • S CT ,/), are derived from fourth-order 
perturbation expansions in the strong coupling limit of the Hubbard model [26]. It is well 
known that such many body exchange interactions play important role in Wigner crystal 
and 3 He solid. Recently, low dimensional quantum spin systems with these terms are well 
studied. In \.J rung > J r i ng = J\ eg case, this model can be solved exactly using matrix product 
method [23]. On this region the ground state is a product of rung-singlet [23] [27]. 

III. OVERVIEW OF EDGE STATES 

It is well known that in finite size Haldane systems, under the open boundary condition, 
quasi four- hold degenerate states appear. This was pointed out by Kennedy [28], who 
calculated the ground state and some lower excitation energies numerically. 

In the Haldane system, such as valence-bond solid (VBS) states[29], a ground state is 
unique in the infinite system. But under the open boundary condition, it is four-hold 
degenerate, which is between the ground state and triplet excitations, and other states, for 
example quintet excitations, have a finite energy gap. 

On the other hand, in the finite system, under the open boundary condition, it is quasi 
four-hold degenerate, which consists of the ground state and triplet excitations. This appears 
universally in Haldane gap systems. These quasi energy gaps between four states decay 
exponentially, as increasing system sizes. This is called as edge states. This can be confirmed 
experimentally [30]. Edge states should appear in not only S = 1 chain but also S — 1/2 
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two-leg ladder system. 

Usually we can confirm edge states investigating the ground state energies under the 
open boundary condition with higher magnetizations, for example a quintet. However, in 
the finite system, the energy gap between the singlet and triplet remains, which results from 
the finite size effect. 

Therefore, when the singlet-quintet gap is small, it is difficult to investigate the edge 
states directly. 



In this section we present our numerical results for low-lying levels of the model for various 
couplings and system sizes, which are following as, 

• singlet: total S z — ; parity=even 

• triplet: total S z — 0, ±1 ; parity=odd 

• quintet: total S z = 0,±1,±2 ; parity=even. 

Here we consider the parity about the inversion symmetry for the leg direction. In the 
commensurate region, the singlet is the ground state. 

Now we study J\ eg = J rung = 1 case with verifying J r i ng - We consider spin gaps between 
the singlet and the triplet and between the singlet and the quintet. 



We calculate spin gaps, using a finite size DMRG algorithm with open boundary condi- 
tion. We use the maximum system size L — 112 and m = 300 in this paper, where m is the 
number of the bases for the truncated left-half system. 



Here we consider a ratio of the energy gap of the triplet-singlet and the energy gap of the 
quintet-singlet, 



IV. NUMERICAL RESULTS ON PHASE TRANSITION 



AE st = E (triplet) — E (singlet) , 
AE sq = E (quintet) — E (singlet) . 



(3) 
(4) 



A. 



RATIO OF ENERGY GAPS 




AE. 
~AE 



sq 



st 



(5) 
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In a normal rung-singlet phase (see Fig.3), this ratio is about 2, based on standard 
magnon picture. If the interaction of magnons is repulsive, then this value should be 2 in 
the infinite system, and it approach 2 from a large value increasing the system size in the 
finite system. In this phase, the triplet excitation state has one magnon, and the quintet 
excitation state has two magnons. If two magnons are independent, this ratio is equal 2. 
This is true in enough large systems (see Fig. 5). 

In a massless phase in a finite system under open boundary condition, excited energies 
are given as 

AE(L) = ^x s (6) 

where x s is a surface critical exponent, v is a Fermi velocity, and L is a system size based 
on CFT. If a phase is massless, then this ratio converge to a finite value x' s /x s . In reality, 
there remains the finite size effect from logarithmic correction(see Appendix) 

In a staggered dimer phase, A st — > and A sq — > finite, so / — > oo. 

So we can distinguish a massless phase from some degenerate phases in large J r i ng region, 
as follows, 

• Massless phase : / — * converge (finite) 

• Staggered dimer phase : / — > diverge 

We find that this ratio tend to diverge in large J r i ng region (see Fig. 5, 6). This shows that 
/ — > 2 in the small J r ing region, and / — > oo in the large J r ing region. So we can deny that 
there exists an extended massless region in the large J r ing region. And it shows that the 
interaction of magnons is repulsive in the rung-singlet phase. We can understand / — > 2 to 
consider that the rung-singlet phase appears. And we can interpret / — > oo as A st — > and 
A sq — > finite. Considering SU (2) symmetry, we can conclude that the edge state appears. 

Here we consider the function /. This phase transition is second order, where the critical 
point is massless, so we can expect that this function behave 

• Jring < Jring ■ f ->fillite ~ 2 

• Jring = J ring '■ f ~> X s (sq) / X s (st) ~ finite 

• Jring > J ring '■ f ~ > 00 
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as increasing system size L, where x s (st), x s (sq) are scaling dimensions. This function is 
independent of system size L at J r i ng = Jring- Therefore, when we plot this function for 
J ring i we can expect that this function behaves as the traditional finite size scaling function 
LAE. This is based on a phenomenological renomalization group analysis [31]. 
In Fig. 7 we define J^ng 3 (L) as 

/ (L + 20, JZT (L + 10)) = / (L, JZT (L + 10)) (7) 
Then, we extrapolate J^™ 3 as follows, 

JZ7 (L) = J™ (oo) + a± + b± + higher order (8) 
Here we neglect higher order terms, so we obtain J%™ 8 (°°) — 0.19379 (see Fig.8). 

B. THE CORRELATION LENGTH 

In this section, we discuss the correlation length in the rung-singlet phase. We determine 
the correlation length £ according to the following formula, 

/ (L, J ring ) = / (L = OO, J r i ng ) + aK (L/f) (9) 

where a is an amplitude, is the modified Bessel function, and / (L — oo, J r j ns ) is about 2. 
As is well known, the dimensional quantum system can be mapped onto d+ 1 dimensional 
classical system[32]. If the one-dimensional quantum system has the energy gap, it is called 
as massive, and its correlation function becomes of the 1 + 1 Ornstein-Zernicke form, that 
is the modified Bessel function, 



/ 



exp[iq-l] d J exp[-|l|/f] : d= 1 
2 , c 2 d q oc I 10 
« +^ ^o(|l|/0 :^ = 2 . 



This form is supported by Quantum Monte Carlo simulation [33] and DMRG calculation 
[34] 

Using the nonlinear least squares method, we determine f(L = oo, Jring)-, a, £• We obtain 
l/£ as a function of Jring-, fitting for L = 20 ~ 112. This figure shows clearly that the l/£ is 
proportionate to J ring- Note that in the case of J ring = 0.175, £ is larger than the system size 
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L. The correlation length diverges at the second order phase transition point. Our result is 
consistent with this fact. 

Near the critical point described as k = 2 SU(2) WZW model, l/£ behave 

— OC | Jring ~ J r ing\ > (H) 

where J£ ing is J ring at the critical point. Using the method of least squares, we obtain the 
£ _1 = at J^ing = 0.1943 (see Fig. 9). This value is close to previous results [16, 19]. 

C. SCALING DIMENSION 

In previous section, we discussed the divergence of the correlation length. In this section, 
returning to the CFT analysis with the periodic boundary condition, we discuss a scaling 
dimension independently. In our previous paper [16], we obtained a scaling dimension x — | 
thus we have calculated k = 2 WZW type central charge c = |. But this model, level k = 2 
577(2) WZW model, has another scaling dimension x — 1. This scaling dimension is related 
with the energy gaps (the correlation length oc AE) near the critical point, 

AE oc \j rmg - J c rmg \ V (12) 

where rj — 1. Miiller et al [18] got this critical exponent, using the cluster expansion and Pade 
approximation. This is based on rung-dimer limit J rU ng oo. However we are interested in 
Jieg = Jmng case. Now we consider a scaling dimension based on CFT. The relation between 
an excitation energy and a scaling dimension , under the periodic boundary condition, is 

27T7V 

Ei-E = — Xi (13) 

where L is a systems size, v is a Fermi velocity in the system and Xi is a scaling dimension. 
In addition to the ground state energy and the excitation energy, we need to obtain v 
numerically. This velocity is obtained from 



2tt 
1 



E[q = — )-E gs 



(14) 



where q is a wave number. Then we extrapolate v(L) as follows, 



1 . / 1 ^ 2 



v (L) =v 00 + a— + b I — ) + higher order (15) 
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We neglect higher order terms and use v^. 

Unfortunately there are logarithmic corrections from the marginal operator in k = 2 
577(2) WZW model. So we must remove corrections, when we calculate scaling dimensions. 
After removing them with the method in Appendix, a scaling dimension is almost inde- 
pendent on the system size at the critical point J r ing = 0.192 (see Fig.10). This result is 
consistent with previous results [16, 18]. 

V. CONCLUSIONS 

In this paper we have investigated the ordered state of a S — \ quantum spin two-leg 
ladder system with the cyclic four spin exchange using DMRG method in the open boundary 
condition and an exact diagonalization method in the periodic boundary condition. 

Generally speaking, it is very difficult to determine whether some value is zero or not 
numerically, especially near the critical point, where the correlation length becomes large 
compared with the system size. 

Considering the ratio of energy gaps of the singlet-triplet and the singlet-quintet, we find 
that the singlet-triplet is degenerate and the singlet-quintet is not degenerate after the phase 
transition. This result supports that the staggered dimer phase appears. In the staggered 
dimer phase with the open boundary condition, since the edge state appears, the singlet- 
triplet energy gap decays exponentially. On the other hand, the singlet-quintet gap is finite 
in this phase. 

We have also investigated the correlation length near the critical point. That result is 
consistent with previous results [16, 18, 19]. 

Honda et al [20] concluded that the massless phase appeared, using DMRG with the open 
boundary condition. But we can consider that they found the small edge state. We obtain 
the unified interpretation which is consistent with previous results [16, 18-20, 22]. 

This problem is related to electron systems on a two-leg ladder. Our result is consistent 
with Tsuchiizu et al [35]. They discussed the ground state phase diagram of half filled 
two-leg Hubbard ladder with inter-site Coulomb repulsions and exchange coupling. They 
insisted that the universality class of the phase transition between the rung-singlet phase 
and the staggered dimer phase is c = 3/2 CFT or the first order phase transition. Our result 
is consistent with their result. 
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According to Hikihara et al this model has a self-duality (J r i ng = J[ eg /2), which is between 
a spin and a vector chirality [22]. In the case that the system has duality and there exists 
a phase transition, if a phase transition point is not the self-dual point, there must be 
another phase transition. Now Hamiltonian is not invariant for the duality transformation 
exactly. But we can expect that there is another phase transition reflected in the duality 
transformation. This should be a future problem. 
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APPENDIX A: WESS-ZUMINO-WITTEN MODEL AND LOGARITHMIC 

CORRECTION [17, 36, 37] 

The minimal conformal field theory with the smallest spectrum containing currents obey- 
ing the Kac-Moody algebra with central charge c is the Wess-Zumino-Witten non-linear a 
model, with topological coupling constant k. Thus the relation between c and k is 

k dim G , . _ s 

C ^JTW (A1) 

where dim G is dimension for an arbitrary representation and -± is called as the dual Coxeter 
number. In the present case, since the system has SU(2) symmetry, we have dimG = 3, and 

Ca o 

2 ~ Z - 

Now we think of the case with the SU(2) symmetry. In this case, primary fields can be 
classified according to their left and right moving spin. There are operators with sl = sr = 
0, h, §. Their scaling dimensions are 

(A2> 

The operators with sl half-odd integer (integer) are odd (even) under translation by one 
site; thus they correspond to states with momentum n (zero). 
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Here we think of the case of the k — 2 SU(2) WZW model. In the case of s L — s R — |, 
which forms \ ® \ = © 1, the scaling dimension is x — §. In the case of s L = s R = 1 
which forms 1©1 = 0©1©2, the scaling dimension is x = 1. This is related to the critical 
exponent of a energy gap. 

Besides Kac-Moody primary operators, there is a marginal operator for all k, namely 
Jl • Jr (which is a primary field with respect to the Virasoro algebra). Current operator 
J L has the conformal weight (h, h) = (1, 0) and J R has the conformal weight (h, h) = (0, 1). 
Thus the J L ■ J R operator has the conformal weight (h, h) = (1, 1) which has the scaling 
dimension x — h + h — 2 and the conformal spin 0, i.e. wave number q — 0. 

Current operators Jl, Jr themselves have conformal spin ±1, corresponding to the wave 
number q = ±2n/L, thus they are related with the spin wave velocity. 

According to the non-Abelian bosonization [38], in SU(2) symmetric gap less system, the 
excitation energies (for q = states) in the periodic boundary condition, including log 
corrections, are 

at? t? t? 2irv ( (S L -S R )\ 

AE = Ei - E w — ( Xi - lnL J ■ ( A3 ) 

where L is a system size, and v is a spin wave velocity. Here S = + S R is total spin which 
is a conserved quantity. Note that 

Sl • S R = \ (S L + S R ) 2 - is| - (A4) 
= + 1) - \sl (s L + 1) - \s R (s R + 1) . (A5) 

In the case of sl — s R — 1, since sl © s R —l © 1=0 © 1 © 2, it forms a quintet s = 2 triplet 
s = 1 and a singlet s = 0. The states sl = s R = 1 correspond to a wave vector q = 0. 
First s = (singlet) case, 

{singlet \ Sl • S R \ singlet) = —2 (A6) 

Second s = 1 (triplet) case, 

(triplet | Sl • S R \ triplet) = — 1 (A7) 

Third s = 2(quintet) case, 

(quintet \ Sl • S R | quintet) = 1 (A8) 



Therefore we obtain 



and 



and 
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AE(s = 2) = Ei - E « ^ ( Xi - ^j- 1 , (A9) 



27Tf / 1 

AE(s = 1) = ^ - E « — f + — ) , (A10) 



27Tf / 1 

AE( S = 0) = E,-Eo^ — (zi + 2— J , ( A11 ) 



thus we can remove logarithmic corrections of energy gap. 

AE(s = 2) + AE(s = 1) = 2—- Xi , 

Lj 



or 



Xi = [AE(s = 2) + AE(s = 1)] . (A12) 
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FIG. 1: Two-leg spin ladder with four spin exchange. 
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FIG. 2: Four spin cyclic exchange on the two- leg ladder. 
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FIG. 4: A staggered dimer order state which consists of two spins S Qj i, S Qj i+i enclosed by a dotted line, 
the open boundary condition. This phase has an order parameter {S a ,i ■ S a ,i+i) -(S a ,» • S Q ,i+i) 



FIG. 3: A rung singlet state which consists of two spins Si^, enclosed by a dotted line. 
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FIG. 5: Ratios of energy gaps as a function of the system size L 




FIG. 6: Ratios of energy gaps as a function of the system size l/L 
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FIG. 8: J^l° s g s as function of the system size 1/L. 
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FIG. 9: The correlaiton length as a function of J r i ng is computed using the modified Bessel function. The 
dot line was computed using a least-squares method. 
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FIG. 10: Size dependence of the scaling dimension x = 1 after removing logarithmic correction. 



